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SUMMARY 


f This report presents an evaluation of the transonic— wing— analysis computer code 

TWING. TWING utilizes a fully implicit approximate factorization iteration scheme to 
solve the full potential equation in conservative form. A numerical elliptic-solver 
grid-generation scheme is used to generate the required finite-difference mesh. 

Several wing configurations have been analyzed, and the limits of applicability of 
the code have been evaluated. Comparisons of computed results have been made with 
available experimental data. Results indicate that the code is robust, accurate (when 
significant viscous effects are not present), and efficient. TWING generally pro- 
duces solutions an order of magnitude faster than other conservative full potential 
codes using successive— line overrelaxation. The present method is applicable to a 
wide range of isolated wing configurations including high-aspect-ratio transport wings 
and low-aspect-ratio, high-sweep, fighter configurations. 


I . INTRODUCTION 


Recently, computers have had an increasing role in the aerodynamic design of 
aircraft for transonic flight. The computation of three-dimensional transonic flows 
about wings in many senses has become a routine affair. However, currently available 
methods still are unable to accurately and efficiently predict transonic character- 
istics of wings with large variations in aspect ratio, Mach number, angle of attack, 
sweep, taper ratio, etc. 

During recent years, inviscid transonic flow computations have been obtained 
using two formulations: (1) the transonic small disturbance (TSD) formulation and 

(2) the full potential (FP) formulation. The TSD formulation is valid for isentropic, 
irrotational flows involving thin bodies at small angles of attack and for transonic 
conditions. Flow-tangency boundary conditions are imposed on mean geometrical sur- 
faces, for example, wing "slits" or fuselage "boxes." Thus, simplified grids can be 
used, usually involving sheared, stretched Cartesian coordinates. This is the primary 
reason why the TSD formulation was introduced before the FP formulation and why more 
complete geometrical configurations have been solved via the TSD approach. 

The worst problem associated with the TSD formulation occurs at the wing leading 
edge where the stagnation region generally causes the small disturbance assumption to 
break down. However, this formulation yields results which are generally better than 
expected. Indeed, good correlations with experiment have been obtained for a wide 
variety of applications. Example applications using the TSD formulation can be found 
in references 1-5. Although the TSD formulation is capable of computing flows for 
relatively sophisticated geometries , the FP formulation is generally desired when an 
accurate solution near the wing leading edge is required. 
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The FP formulation is valid for isentropic, irrotational flows ranging from 
incompressible to transonic conditions. The flow-tangency boundary condition is pre- 
scribed in an exact manner (relative to the TSD approach) at the geometry surface. 

Thus a more complex mapping procedure is required with this approach. Example appli- 
cations using the FP formulation can be found in references 6-10. In addition, a 
survey of computational transonic methods in which the TSD and FP formulations are 
compared and contrasted is presented in reference 11. 

Recently, a large number of researchers have investigated various ways of improv- 
ing solution convergence. This research direction has been largely necessitated by 
the slow rate of convergence associated with the standard relaxation algorithm, 
successive-line overrelaxation (SLOR) . Some of the new algorithms investigated 
include multigrid (refs. 12-16), the strongly implicit procedure (SIP) (ref. 17), and 
various approximate factorization schemes (refs. 18-23). 

In this study, the TWING computer code (transonic wing analysis) has been used 
to calculate the flow field around a variety of isolated wing geometries. The intent 
of this report is to provide a quantitative evaluation of this code with respect to 
reliability and accuracy in solving transonic-wing flow fields. This code uses the 
AF2 fully implicit approximate factorization scheme to solve the FP equation in con- 
servative form (ref. 24). Supersonic regions of flow are stabilized in the spatial- 
differencing scheme by using an upwind bias in the density evaluation (refs. 19, 25, 
and 26). This causes the spatial differencing scheme to be first-order accurate in 
supersonic regions and second-order accurate in subsonic regions. 

Variations in wing aspect ratio, taper ratio, sweep, twist, and airfoil cross 
section have all been investigated in the present study. General details of the TWING 
solution algorithm are discussed first, followed by detailed grid and flow-solver 
results. Correlations with experimental results are included for most of the cases 
presented. 


II . NUMERICAL ALGORITHM 


Governing Equations 

The three-dimensional full-potential equation written in strong conservation-law 
form is given by 


P = ~ ^ T~ ? ' + <!>^ + 4>2)1 

L Y + 1 X ^y ^z J 

The density p and velocity components <py, and (p^f are nondimensionalized by the 

stagnation density pg and the critical sound speed a^, respectively; x, y, and z 
are Cartesian coordinates in the streamwise, spanwise, and vertical directions, 
respectively, and y is the ratio of specific heats. 

equation (i) is transformed from the physical domain (Cartesian coordinates) to 
the computational domain by using a general independent variable transformation. 

This transformation (see fig. 1), indicated by 


(la) 

(lb) 
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C = 5(x,y,z) 
n = n(x,y,z) 
C = C(x,y,z) 


( 2 ) 


maintains the strong conservation-law form of equation (1) (see ref. 27). The full- 
• potential equation written in the computational domain (5 - n - ? coordinate system) 
is given by 
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Note that this formulation for the FP equation is more general than the formulation 
24 . The simplification that all r) = constant surfaces coincide with 
y = constant planes has been removed. Thus, completely general orientations of the 
C - n - 4 mesh can be supported with the present formulation. For more discussion 
on the details of this transformation procedure see reference 24 . 


Grid Generation 

The grid generation scheme used in the present three-dimensional formulation is 
a simple extension of the two-dimensional scheme presented in reference 28 . The • 

finite-difference mesh is generated using a standard two-dimensional algorithm. This 
requires solution of two Laplace equations in each spanwlse plane used as a defining 
station. These equations are transformed to the computational domain and solved 
using a fast approximate-factorization relaxation algorithm (ref. 28 ). 

This establishes values for x and z in each spanwlse plane used as a defining 
station. Coordinate values (x and z) for computational planes between the defining 
stations are obtained via linear interpolation. For the case of a wing with no 
taper or section variation, only two defining stations are required, one at the root 
and one outboard of the tip in the wing-extension region. The root station is user- 
specified, but the wing-extension station is always chosen as a flat plate. In addi- 
wing taper, twist, thickness, and sweep variations can be specified at each 
defining station. 

The coordinate values in the spanwlse direction (y values) are computed from a 
stretching formula, that in its simplest form gives equal spacing over the wing with 
relatively rapid stretching beyond the tip. The x and z values generated for the 
first station outboard of the tip are used for each wing extension station. Then a 
smoothing step for the two y = constant planes on either side of the tip is imple- 
mented to remove grid discontinuities, which may arise in the transition from an air- 
foil of positive thickness to one of zero thickness. An example grid generated using 
the procedure just discussed will be presented in the section on computed results. 


Spatial Differencing 

A finite-difference approximation for equation ( 3 ) suitable for both subsonic 
and supersonic flow regions is given by 

where the operators ), ), andl^( ) are first-order-accurate, backward- 

operators^in the h» and ? directions, respectively. The density 

coefficients p and p are defined by 


i+i/2,j,k ^^*^^1+1/2, j ,k ^1+1/2, j ,k'^i+r+i/2,j ,k 

i,j+i/2,k “ ~ '*^^^i,j+i/2,k '^i,j+i/2,k'^i,j+s+i/2,k J 


( 7 ) 
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where 
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The quantity Mi,j,k is the local Mach number, C is a user-specified constant 
(usually between l!o and 2.0); and the quantities U, V, and W are the contravariant 
velocity components computed using standard, second— order— accurate finite— difference 
formulas. The density p is computed from the second-order-accurate discretized 
version of equation (3b) and is stored at half points in the finite-difference mesh 
(i.e., at i+l/2,j,k). Values needed at i,j+l/2,k or i,j,k+l/2 are obtained using 
simple four-element averages. 


The spatial differencing scheme given by equations (6-9) is centrally differ- 
enced and second-order accurate in subsonic regions. In supersonic regions, the dif- 
ferencing is a combination of (1) the second-order-accurate central differencing used 
in subsonic regions and (2) the first-order-accurate upwind differencing resulting 
from the upwind evaluation of the density coefficients (p and p). In the present 
formulation upwind evaluations of the density have been used along only the i and j 
indices (wraparound and spanwise directions, respectively). Thus the fully rotated 
form of the spatial difference scheme has not been used. This has not adversely 
affected convergence for any of the cases presented herein. However, other cases 
involving large regions of supersonic flow at the wing trailing edge may require the 
fully rotated difference scheme with all three values of density in equation (6) 
biased in the upwind direction. 


AF2 Iteration Scheme 


The AF2 fully implicit iteration scheme used in the present study can be 
expressed in a three-step format given by 


Step 1: 


Step 2: 


Step 3: 


(a - ^ A,^ |g? . = auL(|)? . , + cxA, f^ . 

\ j ny®i,j ^i.J.k k+i i,j. 


i ,k+x 




(“-*-V^i,3,k = ^?,J,k 


(10a) 


(10b) 


(10c) 
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where the n superscript is an iteration index; a is an acceleration parameter; 
^^i,j,k iteration residual (defined by eq. (6)); o) is a relaxation 

factor, equal to 1.8 for all cases presented; gi^j is an intermediate result stored 

at each grid point in a given k plane; and fi,j,k intermediate result stored 

at each point in the finite-difference mesh. The Ai, A j , and A^, coefficients are 
defined by 
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(pAj/J) 
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i,j ,k-i/2 




( 11 ) 
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and the density coefficients p and p are defined by equation (7). 


The quantity a appearing in equation (10) can be considered as At“^. The best 
choice for a is to use a sequence of values. The small values are particularly 
effective for reducing the low-frequency errors, and the large values are particularly 
effective for reducing the high-frequency errors. The a sequence given in refer- 
ence 19 has been used for the computations presented herein. 


A ({i^t-type term has been added inside the brackets of step 2 (see eq. (10b)). 
This term is necessary to provide time-like damping to the iteration scheme. The 
parameter 3^ is user-specified and usually ranges between 0.05 and 0.3. Larger 
values are generally required for cases with larger regions of supersonic flow. 


In step 1, the g array is obtained by solving a tridiagonal matrix equation for 
each 5 = constant line in the kth plane. In step 2, the f array is obtained from 
g by solving a tridiagonal matrix equation for each n = constant line, again for 
just the kth plane. Next, step 1 is used to obtain the g array for the k - 1 
plane, and then step 2 is used to obtain the f array for the k — 1 plane, etc. 

This process continues until all values of f in the three-dimensional mesh are 
established. Then, by using step 3, the correction array is obtained from the f 
array by solving a simple bidlagonal matrix equation for each C line in the entire 
finite-difference mesh. Other aspects. of the present spatial differencing and itera- 
tion schemes are discussed in more detail in reference 24. 


, , Bbund^^ Conditions 

’ -:The wing- 5 urf ace houndhry? Conditi^^ body-fitted coordinate system is 
implemented .by.vapplyirig 't *' -fr ■ 


^ _ /pw\ 

V? / i , j ,NK+i /2 \ /i , j ,NK-l 


/2 


( 12 ) 


where NK is the wing, surface. A similar reflection condition is used for the sym- 
metry plane (y = 0). Along the free stream and outer boundaries the initial free- 
stream distribution of ' (j) is fixed for nonlifting conditions. For lifting wings the 
outer boundary is updated by the usual compressible vortex solution with circulation 
Fj. A user-specified relaxation parameter (RGAM) is used to control the circulation 
build up. Details describing the circulation algorithm can be found in reference 24. 
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III. COMPUTED RESULTS 


In this section several results computed from the TWING computer code are pre- 
sented. In each case the solution procedure started with the numerical generation of 
the finite-difference mesh. An example grid is shown in figures 2-9 for the Wing C 
geometry obtained from reference 29 (configuration 4 in table 1). The x/c, y/c, 

Ld z/c Cartesian coordinates used in these figures are based on an oripn at the 
wing-root midchord and are normalized by root chord c. This grid contains 
127x27x20 = 68,580 points (wraparound, spanwise, and near-normal directions, respec 
tively), and is typical of the grids used for the results presented herein. Note 
that with this grid 127x17 = 2159 points are used to surface geometry. 

This configuration is a high-sweep (Ale = ^5°), high-taper (TR 0.3) wing wit 
twist. The amount of taper and twist is apparent from figures 2, 3, and 4 which show 
blowups of the root, tip, and wing extension airfoil sections plotted to the same 
scale. Note the efficiency with which the "0 mesh" topology clusters points around 
the airfoil-section geometry regardless of its position in the span (root or tip). 

For this particular grid a special clustering procedure has been used to cluster grid 
points at the wing leading edges. This is especially useful for wings with sharp 
leading edges. Figures 5 and 6 show additional detail about the wing root section 
leading and trailing edges, respectively. The wing planform including the surrounding 
mesh is shown in figure 7, and the grid in the 5 - constant ° 

the wing half-chord position (i.e., the wing on edge) is shown in figure 8. Apparent 
in these views is the wing-tip grid topology and the stretching rates used above and 
below the wing as well as outboard of the wing tip. Finally, figure 9 contains a 
perspective view of the three-dimensional grid which consists of the symmetry plane, 
wing surface and vortex-sheet grid distributions. This view very much represents the 
view of the wing mounted in the wind tunnel (except for the vortex sheet) and presents 
a good picture of the overall grid topology. 

The first case presented consists of a large aspect ratio wing (AR =40) with 
zero sweep. The airfoil section (Korn airfoil, see ref. 30) is constant in the span 
wise direction. The freestream Mach number is 0.75 and the angle of attack is 0.115 . 
This case is basically two dimensional in nature (everywhere except at the tip) and, 
thus, is compared with the two-dimensional hodograph result of reference 30 in 
ure 10. The agreement for this transonic "shock free case is quite good. The TWING 
code produces a small shock at the aft sonic line which is similar to other results 
produced for this sensitive problem (ref. 31). 

The next case involves a subcritical calculation (M„ = 0.7, a =4°) about a 
moderately swept, tapered wing without twist (see conf. 1 in table 1). A comparison 
of wing surface pressures with experiment (ref. 32) is shown in figure 11. The com 
putatiLal angle of attack was set so as to approximately match the experimental value 
of lift rather than the experimental value of the angle of attack. This is primari y 
required to simulate viscous effects which are not modeled in the inviscid TWING 
formulation. Overall the agreement between the experiment and the computational 

result is excellent. 

The next case involves the "Wing A" geometry of reference 29 (see also ref. 33). 
Characteristics of this high-taper, low-sweep wing are presented in table 1 (conf. ). 
A comparison of wing surface pressures with experiment is shown in figure tor 
measured experimental conditions of M„ = 0.818, a = 2.9°. Corrections to the experi- 
mental Mach number and angle of attack of -0.02 and -2.3° were required in the present 

computation (i.e.. M„,com = 0-798, and a^om = 0.6 ) to produce good correlation with 
the experimental pressure distributions. The wind tunnel wall-interference study 
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conducted in reference 29 produced similar, but smaller corrections for this set of 
conditions (AM = -0.005 and Aa = -1.4°). The angle of attack correction of refer- 
ence 29 was determined from the FLO 22 computer code (ref. 34) by varying the compu- 
tational angle of attack until the computed upper and lower surface pressures in the 
wing leading edge region matched experiment. Then the difference between the experi- 
mental and computational angles of attack was taken as the desired angle of attack 
correction Aa. After Aa was obtained, the Mach number correction (AM) was deter- 
mined by using the Bailey-Ballhaus TSD computer code (ref. 35). Experimentally mea- 
sured wind-tunnel-wall pressures were imposed as boundary conditions on the outer 
boundary of the computational domain (set to model the wind-tunnel wall position) . 
Next, free-air solutions were computed with variations in the Mach number to produce 
a "best match" free-air solution with the wind tunnel solution. The difference in 
Mach number (if sufficiently small) between the wind tunnel wall case and the "best 
match" free air case was taken to be AM. 

This wind tunnel correction procedure seems to be a good one as it includes 
nonlinear, three-dimensional effects, but could suffer from numerical uncertainties 
arising from viscous modeling as well as other numerical errors, for example, noncon- 
servative differencing, mesh effects, etc. It should be noted that the corrections 
determined in reference 29 for the Wing A geometry, as well as other geometries to be 
discussed shortly, produced good agreement for comparisons between FL022 and experi- 
ment, but not good agreement between other computer codes and experiment. The 
probable reason is that this correction procedure removes two types of errors: 

(1) wind-tunnel wall Interference errors and (2) numerical errors associated with each 
computer code. This suggests that separate corrections tailored to each Individual 
computer code may be required. 

In the present case the Increased angle of attack correction (Aa = -2.3°) is 
probably the result of viscous effects but the discrepancy between the Mach number 
corrections is puzzling. At any rate, with the appropriate corrections, the present 
results are in good agreement everywhere except near the base of each shock outboard 
of the wing root where obvious shock/boundary layer Interactions affect the local 
pressure distribution. The last pressure distribution at y/c = 0.95 is in only 
fair agreement with experiment and probably suffers from numerical modeling inade- 
quacies at the tip. 

The moderate sweep case presented in reference 29 (Wing B) is discussed next. 

See table 1 for specific characteristics of this geometry (conf. 3). Pressure coeffi- 
cient comparisons with experiment are presented for Wing B at M„ = 0.9, and 
a = 3.9° (experimental conditions) in figure 13. The AM and Aa corrections used in 
the TWING results were —0.02 and —1.7°. The corresponding AM and Aa corrections 
computed in reference 29 were -0.005 and -1.0°. Overall, the agreement between 
experiment and TWING is good at every span station. The aft shock position near the 
root is about 5-10% too far downstream, but moves into good agreement outboard of 
the midchord position. 

The last configuration presented from reference 29 is a high-sweep, low-aspect- 
ratio wing — Wing C. The characteristics of this geometry are presented in table 1 
(conf. 4) and have already been displayed in figures 2-9 where several views of a 
sample grid about this configuration are presented. Three flow conditions for Wing C 
are considered in the present study: a weakly supercritical case, M^ = 0,7, 

a = 4.949°; a moderately supercritical case, M„, = 0.85, a = 5.9° (design condition); 
and a strongly supercritical case M„ = 0.9, a = 4.914°. Pressure distribution com- 
parisons between the present computed results and experiment are shown for the first 
case in figure 14. For this case the experimental Mach number was matched and the 
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angle of attack correction was -0.55®. Except for a slight overprediction of the 
minimum pressure at the leading-edge upper surface the overall agreement 'is excellent. 

The second Wing C result is compared with experiment in figure 15. The Mach 
number and angle of attack corrections for this case are AM = -0.02 and Aa = -0.9®. 
The agreement for this case is quite good everywhere except at the tip where the need 
of viscous corrections is apparent. Of particular note in this calculation is the 
ability of the TWING code to predict the oblique shock which exists at both the third 
(y/c = 0.5) and fourth (y/c = 0.7) span stations. The differencing scheme in this 
region is entirely first-order accurate and yet little shock smearing is exhibited. 

The final Wing C calculation (M^, = 0.9, a = 4.914°) is compared with experiment 
in figure 16. For this case the numerically determined Mach number correction of 
reference 29 was AM = -0.01, and the corresponding AM used in TWING was -0.02. 
Overall, the agreement between these results is good. Viscous effects probably cause 
the discrepancy in shock strength and position at the tip station (y/c = 0.9). 

The experimental and computational spanwise load distributions for the three 
Wing C cases just presented are compared in figure 17. The agreement is best for the 
lower Mach number case and worsens somewhat for higher Mach numbers. This is probably 
due to the larger angle-of-attack corrections required for the higher Mach number 
cases to match the pressure distributions in a reasonable way. Nevertheless, the 
agreement is good in each case. 

A convergence history for the last Wing C calculation is shown in figure 18. 

Both the lift and the number of supersonic points (NSP) plotted versus iteration num- 
ber n, and computer time in seconds on the CRAY-IS computer are displayed. The lift 
and NSP levels are essentially established in about 20 iterations. After about 
40 iterations the lift and NSP levels stabilize and the solution is essentially con- 
verged. More interesting is the amount of computer time required to establish the 
solution. The entire convergence history shown in figure 18 required only about 
25 sec of computer time and the crudely converged solution (40 iterations) required 
about 13 sec. Solutions obtained with codes using the SLOR algorithm, similar 
meshes, and the same computer would require at least an order of magnitude more com- 
puter time to achieve a similar level of convergence. 

The small amount of computer time required per iteration, which is just as impor- 
tant as the small number of iterations required for a solution, is caused by several 
factors. The first is the high speed of the CRAY-IS computer. Second, the small 
number of operations packed into the fully implicit AF2 iteration scheme represents 
a minimum burden to any computer. Last, and most important, the way in which the 
operations of the AF2 iteration scheme lend themselves to efficient vector processing 
greatly speeds the computational processing rate within the TWING computer code. 

Other algorithms, such as SLOR, would not be so efficiently vectorized because of 
inherently recursive operations. 

To obtain an idea of the effects of the angle of attack (Aa) and Mach nvunber (AM) 
corrections discussed earlier, several pertinent results are shown in figure 19. The 
pressure coefficient distributions at two span stations for configuration 4 (n ~ 0.30 
and 0.70) are compared with experiment for several conditions: (1) the uncorrected 

experimental conditions (Moo = 0.85, a = 5.9°), (2) the corrected experimental condi- 
tions using the corrections cited in reference 29 (M^ = 0.845, a = 5.0°), and (3) the 
corrected experimental conditions using the present corrections (M^„ = 0.83, a = 5.0°). 
As seen from figure 19 the angle of attack correction is more important than the Mach 
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number correction. The conditions cited in reference 29 yield a reasonable solution 
in the present case, primarily because both angle-of-attack corrections are the same. 

The next two numerical/experimental correlations are for the ATFI (Advanced 
Technology Fighter Integration) configuration obtained from General Dynamics~Fort 
Worth (R, Cox, private communication, 1982), Details of this configuration are given 
in table 1 (conf. 5). This geometry is a high-sweep, low-aspect-ratio wing with an 
unusually small taper ratio (TR = 0.177), Such a TR places heavy demands on the 
numerical scheme because of the highly sheared mesh which results. Two flow condi- 
tions are considered for this configuration: (l)Moo = 0,6, a = 3.99° and 

(2) = 0.9, a = 4,01°. Pressure coefficient comparisons for the first case are 

presented in figure 20, Excellent agreement with experiment is obtained for this 
subcritical case everywhere except at the tip (y/c = 0,897) where the tip missile 
used in the experiment was not modeled in the numerical simulation. The sharp- 
leading edge expansion predicted in the numerical result on the upper surface actually 
extends beyond the plotted results shown in figure 20 to a peak value of Cp -2.7. 

It has been truncated to provide a more accurate scale for the rest of the Cp dis- 
tribution shown at each spanwise station. Repeated refinements of the mesh in the 
leading edge region indicate that the peak is numerically correct and indeed, a 
similar peak, although much smaller in extent, is exhibited at several span stations 
in the experiment. It is anticipated that viscous effects and/or pressure tap aver- 
aging effects are responsible for the mismatch between the peak levels. 

The last ATFI solution considered is an extremely difficult case consisting of a 
freestream Mach number of 0.9 and 4° angle of attack. A comparison between numerical 
and experimental results is shown in figure 21. Except for expected discrepancies at 
the shocks, trailing edge, and tip (because of the absence of the experimentally 
mounted tip missile) the agreement is remarkably good. The leading edge expansion/ 
shock and the aft shock details are accurately predicted at each span station away 
from the tip. 

Figures 22 and 23 present the results of an interesting study to determine the 
ability of TWING to handle both positive and negative sweep. The shock sonic line has 
been plotted in each figure on a normalized planform of the wing. This configuration 
consisted of a parallelogram NACA 0012 wing (i.e., TR = 1.0) at ^ = 0.8, a = 3°. 

The aspect ratio was six. Figure 22 shows results for positive, or aft, sweep 
(sweep = 0°, 20°, 30°, 50°, and 60°), and figure 23 shows results for negative, or 
forward, sweep (sweep = 0°, -20°, -30°, -50°, -60°, and -70°). A case involving 
+70° of sweep was also obtained, but because it was entirely subcritical, no result 
at A = 70° appears in figure 22. As expected, the shock position moves forward and 
decreases in strength for increasing sweep (either positive or negative). For all 
cases the shock approaches the symmetry plane wall (2y/b =0) in a normal fashion 
(i.e., in the physical domain), as it must to satisfy inviscid tangent-flow boundary 
conditions. This causes the variation in position of the shock/ symmetry-plane inter- 
section as shown in figures 22 and 23. 

Of particular interest in these calculations is the robustness displayed by the 
TWING code. The TWING grid mapping becomes singular as the sweep angle approaches 
90°. In addition, cross-derivative terms in the transformed FP equation grow as 
the sweep increases. Since these cross-derivative terms are not represented in an 
implicit fashion in the AF2 iteration scheme, instability might be expected to develop 
for even moderate levels of sweep. However, as indicated by the present results very 
large values of sweep (-70° g A g 70°) can be accommodated by the TWING computer code 
in a stable fashion. 
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A summary of computer times obtained from the TWING computer code for a range of 
computed examples is shown in table 2. The convergence for each of these cases is 
generally nonoptimal. The iteration number corresponds to the number of iterations 
required to reduce the initial maximum residual by two orders of magnitude. Computa- 
tion times based on this measure of convergence range from 28 sec to 67 sec. Essen- 
tially converged results are obtained when the lift is bounded by a 2% band about the 
final value of lift. This level of convergence yields computational times ranging 
from 13 sec to 43 sec. The grid generation time is also displayed and ranges from 
1.3 sec (three defining stations) to 3.0 sec (five defining stations). 

These computational times are generally an order of magnitude smaller than other 
conservative full potential codes using SLOR on the same computer (see ref. 24 for a 
comparison between FL028 and TWING). In addition, the advantage of vectorization, as 
provided by the CRAY- IS vector computer, produces over an order of magnitude improve- 
ment in computational speed for TWING relative to the slower CDC 7600 computer (see 
ref. 36 for more discussion on this point). 


IV. CONCLUSIONS 


The TWING three-dimensional transonic full-potential computer program has been 
evaluated using a wide range of test examples. Calculated pressure distributions have 
been compared with experiment for most of the cases presented. The major conclusions 
of this study are as follows: 

(1) Good correlation with experimental pressure distributions have been obtained 
in generally all cases when viscous effects were not important (that is, when strong 
shocks and the associated shock/boundary layer interaction effects are not present) . 

(2) The TWING computer program has proven to be computationally efficient (even 
without optimizing each individual run for minimum run time) and robust (no conver- 
gence problems were encountered for the cases presented herein). 

(3) The TWING computer program has proven to be geometrically general for 
isolated-wing calculations in that a wide range of geometries has been run including 
high-aspect-ratio, low-taper transport-type configurations and low-aspect-ratio, high- 
taper fighter-type configurations. 

(4) Viscous effects need to be Introduced to model shock/boundary-layer inter- 
actions, decambering and displacement thickness effects, and trailing-edge separation 
and wake effects. 
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TABLE 1.- SELECTED WING GEOMETRY PARAMETERS 


Configura- 

Configura- 

Configura- 

Configura- 

Configura- 

tion 1 

tion 2 

tion 3 

tion 4 

tion 5 

(ref. 32) 

(ref. 29) 

(ref. 29) 

(ref. 29) 

a 


Aspect ratio 6.0 8.0 

m 

Taper ratio .5 .4 

(TR) 

Leading edge 37.25 27 

sweep (Ale) 

Root angle of 0 2.76 

attack 
Or - deg 

Tip angle of 0 -2.04 

attack 
6t - deg 

Root thickness- .12 .12 

to-chord ratio 


3.8 2.6 ' 3.8 

.4 .3 .177 

35 45 40 

2.50 2.38 0 

-4.00 -5.79 3.0 

.06 .07 .04 

.06 .11 .04 


(t/c) j 

Tip thickness- .12 .12 

to-chord ratio 

(t/c)t 

^R. Cox, private communication, 1982. 


TABLE 2.- TWING CONVERGENCE TIMES ON THE CRAY-IS COMPUTER (NONOPTIMAL 
CONVERGENCE), ALL TIMES ARE FOR A 127x27x20 (68580) MESH WITH 
127x17 (2159) SURFACE MESH POINTS 


Wing 


M^ 


Ow» 

deg 


CPU time 
grid 

generation, 

sec 


Number of CPU time, 
iterations sec 


CPU time for 98% 
development of lift, 
sec 


Configuration 


2 

0.798 

0.6 

1.3 

164 

51.8 

31.0 

3 

.88 

2,2 

1.3 

170 

54.9 

12.9 

4 

.7 

4.4 

1.4 

90 

27.6 

17.2 

4 

.83 

5.0 

1.4 

99 

30.9 

18.1 

4 

.88 

3.8 

1.4 

83 

28.9 

14.3 

5 

.6 

3.5 

3.0 

102 

31.3 

18.7 

5 

.9 

3.0 

3.0 

196 

66.8 

43.0 
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Figure 2.- Numerically generated finite-difference mesh about configuration 4 
(ref. 29), root station (127x20 grid points). 
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Figure 3.- Numerically generated finite-difference mesh about configuration 4 
(ref. 29), tip station (127x20 grid points). 
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Figure 4.- Numerically generated finite-difference mesh about configuration 4 
(ref. 29), station outboard of the tip (127x20 grid points). 
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Figure 7.- Three-dimensional grid about configuration 4 (ref. 29), planform view, 

27x20 grid points. 
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Figure 8.- Three-dimensional grid about configuration 4 (ref. 29), half chord view, 

27x20 grid points. 


23 



Figure 9.- Perspective view of the grid about configuration 4, 127x27x20 grid points 
(only every fourth point plotted in the wraparound direction). 
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HODOGRAPH SOLUTION (REF. 30) 

TWING SOLUTION 



Figure 10.- Pressure coefficient comparison (airfoil 75-06-12, ref. 30, 

= 0.75, a = 0.115®). 
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Figure 14. ~ Comparison of pressure coefficient distributions, experiment taken from 
Hinson and Burdges (ref. 29), configuration 4, = 45“ , TR = 0.3, JR. = 2.6, 

= 0.7. 
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Figure 16.- Comparison of pressure coefficient distributions, experiment taken from 
Hinson and Burdges (ref. 29), configuration 4, = 45®, TR = 0.3, M. = 2.6. 
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Figure 19.- Concluded. 
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